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Abstract 

The linear eigenvalue problem governing the stability of the mechan- 
ical equilibrium of the fluid in a electrohydrodynamic (EHD) convection 
problem is investigated. The analytical study is one of bifurcation. This 
allows us to regain the expression of the neutral surface in the classical 
case. The method used in the numerical study is a Galerkin type spectral 
method based on polynomials and it provides good results. 
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1 The physical problem 

The investigated physical model is one of the two EHD models of Roberts 
[5], based on the Gross'experiments [B] which are concerned with a layer of 
insulating oil confined between two horizontal conducting planes and heated 
from above and cooled from below. The experimental investigations showed 
that the presence of a vertical electric field of sufficient strength across the 
layer, lead to a tesselated pattern of motions, in a manner similar to that of 
the classical Benard convection [T2] . Gross [BJ suggested that this phenomenon 
may be due to the variation of the dielectric constant of the fluid with the 
temperature. 

In the first model investigated by Roberts [9] the dielectric constant is al- 
lowed to vary with the temperature. The homogeneous insulating fluid is as- 
sumed to be situated in a layer of depth d (the fluid occupies the region between 
the planes z — ±0.5d, which are maintained at uniform but different tempera- 
tures), with vertical, parallel applied gradients of temperature and electrostatic 
potential. The uniform electric field is applied in the z direction. 

The equations governing the EHD convection upon normal mode represen- 
tations are [T^] 

(D 2 -a 2 - s)(D 2 -a 2 - Prs)(D 2 - a 2 ) 2 F = La 4 F - R a a 2 (D 2 - a 2 )F (1) 
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with the boundary conditions on F 
F = D 2 F = D(D 2 -a 2 - Prs)F = 

at z = ±0.5. 

((D 2 - a 2 )(D 2 -a 2 - s)(D 2 -a 2 - Prs) + R a a 2 )(DF ± kaF) = 

(2) 

Here, the unknown function F is the amplitude of the temperature pertur- 
bation 6, i. e. 9 = F(z)e^ lx+m ^ +st \ Pr is the Prandtl number, k = 

^ e 

with e m the value of the dielectric field at the temperature T m — and e the 

electric constant of the solid in z > — , a is the wavenumber, a 2 = I 2 + m 2 , R a 

is the Rayleigh number, L is a parameter measuring the potential difference 
between the planes. 

Roberts [5] investigated only the stationary case, i.e. s = 0, so the eigen- 
value problem consists from an eight-order differential equation 

(D 2 - a 2 fF - La 4 F + R a a 2 {D 2 - a 2 )F = (3) 

and the boundary conditions 

F = D 2 F = D(D 2 ~a 2 )F = ((D 2 — a 2 ) 3 + R a a 2 )(DF ± kaF) = at z = ±0.5. 

(4) 

He found that when the smallest Rayleigh number, i? Qmin = min a R a (a), varies 
from -1000 to 1707.762, L decreases from 3370.077 to 0. 

The second model [S] was also been investigated by Turnbull [T3], [2]. In 
this case, the variation of the dielectric constant is not important but the fluid 
is weakly conducting and its conductivity varies with temperature. 

The eigenvalue equation has the form [12] 

(D 2 — a 2 ) 3 F + R a a 2 F + Ma 2 DF — (5) 

with M a dimensionless parameter measuring the variation of the electrical 
conductivity with temperature. The boundary conditions, written for the case 
of rigid boundaries at constant temperatures, read 

F = D 2 F = D(D 2 - a 2 )F = at z = ±0.5 (6) 

Roberts [9] found that when M is increasing from to 1000, i? amin is increasing 
from 1707.062 to 2065.034. 

Straughan [T^] also investigated these EHD convection problems, develop- 
ing a fully nonlinear energy stability analysis for non-isothermal convection 
problems in a dielectric fluid. 

2 The bifurcation analysis 

The linear stability of the motion or the equilibrium of a fluid in many 
problem from hydrodynamic, electrohydrodynamic or hydromagnetic stabil- 
ity theory is governed by a linear higher-order ordinary differential equation 
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with constant coefficients and homogeneous boundary conditions. The exact 
solution of such equations or, for the case of eigenvalue problems, the exact 
eigenvalue is most of the times impossible to find. That is why, numerical 
methods, usually implying an infinite number of terms, leading however to 
an approximative solution by some specific truncations to a finite number of 
terms, are used. However, the theoretical methods can impose restrictions with 
regards to the numerical results. 

For the considered problem let us introduce the direct method [5] which 
consists in the determination of the eigenfunctions and their introduction into 
the boundary conditions. 

The characteristic equation associated to <j5j> is 



(A 2 -a 2 ) 6 + Ma 2 X + R a a 2 







(7) 



When the characteristic equation has multiple roots the straightforward appli- 
cation of numerical method can lead to false secular points. That is why, these 
cases must be investigated separately. 

Proposition 1. For M = the only secular points are those situated on 

((2n- iyV+a^ 

NS n : R a = ~ , Vn G N. 

a z 

Proof. For M = 0, the characteristic equation ([7]) reduces to 

(A 2 - a 2 ) 3 + a 2 = 0, with a 2 = R a a 2 . (8) 
In this classical case the roots of (fSl) have the form 



Ai, 2 = \/a z + ^-a 2 ei, 2 , A 3 = \J a 2 + ^a^, 
A4 = — Ai, A5 = — A 2 , A6 = —A3 
so the general solution of ([5]) has the form 



F = Ai cosh(Aiz) + Bi sinh(Aiz). 



Replacing the solution F into the boundary conditions © we get the secular 
equation 
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(9) 



with 



tanh(Ai/2), m = A 2 - a 2 , i = 1,2,3. 
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When cosh(Ai/2) 7^ 0, i = 1,2,3, we can rewrite the secular equation as 
A = Ai • A 2 with 



Ai = Ai//imi(A2 - A3) + A 2 Ai2"'-2(A3 - A?) + X 3 fi 3 m 3 (Xl - \%) 

and 

A 2 = AiTOi(A3^3m2-A2At2"i2)+A2m2(Ai^im3-A3^3mi)+A3m 3 (A2^2nT-i-Ai^im2). 

For a > 0, the equations Ai = and A2 = have only null solutions R = 0, 
so no secular points exists on these surfaces. 

The condition cosh(Ai/2) 7^ 0, i = 1,2,3 is not fulfilled only for i = 3, i.e. 
cosh(A 3 /2) = <^ cos(A 3 /2) = <^ A3 = -(2rt - l) 2 n 2 , which implies that the 
secular curve is NS n . And, indeed, the critical values of the Rayleigh number 
R a belong to NS\ only identical to the classical one from Chandrasekhar 
[?]. □ 

The general form of the solution of the two-point problem for the govern- 
ing differential equation is written in terms of the roots of the characteristic 
equation associated with the differential equation. In addition, this form de- 
pends on the multiplicity of the characteristic roots. Introducing the general 
solution into the boundary conditions the secular equation is obtained and it 
depends on the multiplicity of the characteristic roots. As a consequence, the 
secular equation has different forms in different regions of the parameter space. 
Each eigenvalue is a solution of the obtained secular equation, so the eigen- 
value depends on all other physical parameters. The neutral manifolds (the 
most convenient manifolds from the physical point of view), generated by the 
secular equation separate the domain of stability from the domain of instability. 

Let us consider the general case when the roots of the characteristic equation 
Ai, A2, A6 are distinct. Then the general solution of ([5]) has the form 
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F(z) — Ait XiZ . Introducing it into the boundary conditions J6|) we obtain 

i=l 

the secular equation pQ 

A(a,M,R a ) = 0, (10) 

where A is a determinant. Its i-th column has the same form in Ai as any 
other j-th column in Xj . If A^ = Xj , then the i-th and the j-th columns in A 
are identical. Therefore A = 0. In fact, in this situation, (fTUf is not entitled 
to serve as a secular equation and the direct numerical computations will be 
invalid. When M 7^ 0, some particular cases interesting from the bifurcation 
point of view, arise due to the existence of bifurcation sets of the characteristic 
manifold. 

Let us consider the surface So defined by the points (a, M, R a ) — (a, M, a 4 ). 
In this case we have the following result 

Proposition 2. Let us define the surfaces 



(33T3V2T)y/90±10V2T a 3 . 
SilM= 250 ' l=h2 - 
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The surface So H Si, i = 1, 2 is a bifurcation set of the characteristic manifold 
defined by (0). Tie points on Sq Si, i = 1,2 are not secular. 

Proof. If (a, M, i? a ) € So then i? a = a 4 and one of the roots of the charac- 
teristic equation is , for instance, Ai = 0. Assuming that M ^ and a > 0, Ai 
is not a double root of |J7J|. The search of multiple roots reduces then to the 
equation 

A 5 — 3A 4 a 2 + 3A 2 a 4 + Ma 2 = 0. (11) 

No multiple roots of algebraic multiplicity order greater than 2 exists. The 
double roots of (fTTj) must also be roots of its derivative 5z 4 — 9z 2 a 2 +3a 4 = 0. In 

these conditions the possible double roots are A 2 , 3 = A = - — \/90± 10V21a 

only for (a,M,R a ) £ Si, i = 1,2, i.e. the surfaces So H Si, i = 1,2 are 
bifurcation sets for the characteristic equation ©. In the case of multiple 

n 

roots, the general form of the solution of ([5]) is F(z) = ^ Pi(z)e XiZ , where Pi 

i=l 

is an algebraic polynomial of rrii — 1 degree, m, being the algebraic multiplicity 
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of Aj, in our case F(z) = A + (B + Cz)e x z + J2 Ae XtZ . 

2=4 

Formally, the secular equation is deduced from (jTTJJ) , by writing the column 
i for Aj while the columns i + i + 2,..., i + rrii — l are obtained by differentiating 
I, I = 1, 2, rrii — 1 times the i + l-th column of (jTUJ) with respect to Aj + ; and 
then replacing Aj + ; by Aj. 

However, the numerical evaluations show that now secular points exists on 
these surfaces. 

□ 

3 Spectral methods based study 

The second part of our study regards the numerical treatment of the two- 
point problem © - ©■ 

Large classes of eigenvalue problems can be solved numerically using spec- 
tral methods, where, typically, the various unknown fields are expanded upon 
sets of orthogonal polynomials or functions. The convergence of such methods 
is in most cases easy to assure and they are efficient, accurate and fast. Our 
numerical study is performed using a weighted residual (Galerkin type) spectral 
method. 

Introducing the new function U = (D 2 — a 2 )F, the generalized eigenvalue 
problem 

{D 2 ~a 2 ) 2 U = -R a a 2 F-Ma 2 DF, 
(D 2 — a )F = U, (12) 
F = U = DU = at z= ±0.5. 

is obtained. Following ?„ we consider the orthogonal sets of functions 
{0»}i=i,2,...,jv : ^(z) = Jl 0.5 L\(t)dt, verifying 4(±0.5) = 0, 
{•',}, i.i v : Pi{z) = 5 X! Q 5 L*{t)dtds, verifying ft(±0.5) = Dft(±0.5) - 0, 



5 



with L* — Li(2x) the shifted Legendre polynomials on (—0.5,0.5) and Li the 
Legendre polynomials on (—1, 1). The unknown functions from (|12p . U, F, are 
written as truncated series of functions ft, irrespective <fii, i.e. 

N N 
i=l i=l 

The boundary conditions on U and F are automatically satisfied. Replacing 
these expressions in (j!2[) . imposing the condition of orthogonality on the vector 
(ft, <fik) T , k — 1,2, N, we get an algebraic system in the unknown, but not 
all null, coefficients Ui, Fi. The secular equation, written as the determinant of 
the obtained algebraic system, gives us the values of the Rayleigh number as a 
function of the other physical parameters. The smallest values of the Rayleigh 
number for various values of the parameters a and M form the neutral surface 
that separates the domain of stability from the instability domain. All the 
expression of the scalar products resulting in the algebraic system are given in 
[5] for the case of shifted Legendre polynomials on (0, 1), but they are easy to 
adjust to the interval (—0.5. .0.5). The specific choice of basis functions led to 
sparse matrices, with banded sub-matrices of dimension N x N. 

The numerical evaluations of the critical Rayleigh numbers were obtained 
for a small number of terms N (JV = 6) in the truncated series confirming the 
well-known accuracy of spectral methods. We obtained that critical values of 
R a are increasing from 1734.120 to 2082.808 when M is increasing from to 
1000, similar to the ones of Roberts [9]. 

The unknown vector fields from (|12p can also be expanded upon complete 
sequences of functions in L 2 (— 0.5, 0.5) defined by using Chebyshev polynomials 
that satisfy the boundary conditions of the problem. Keeping the above nota- 
tions, the functions (pi, i — 1,2, ..., N are defined by (pi(z) = T*(z) — T* +2 (z) 

and ft, i = 1, 2, ...,N by ft(z) = T* - ^±^T* +2 + ^T* +4 [11] with T*, 

i = 1,2, ...,7V, the shifted Chebyshev polynomials on (—0.5,0.5) defined in a 
similar manner as the shifted Legendre polynomials. All the evaluations of the 
scalar products were based on the orthogonality relation 

Z" - 5 f -c S if % = i 

/ T:{z)T* m {z) W %z)dz=\ f^ ltl h (13) 

•'-0.5 [ 0, if I ft J, 

with respect to the weight function w*(z) = — — 

Vl/4-z 2 ) 

The numerical results where obtained for a larger number of terms in the 
expansion sets (N =11) and they show that the shifted Legendre based method 
is more effective in this case. We can mention that the Chebyshev polynomials 
are considered suitable more likely for the tau method or the collocation type 
methods. Some numerical evaluations of R a as a function of a and M are given 
in Table 1. 

Other sets of complete orthogonal functions based on Chebyshev polyno- 
mials and satisfying various boundary conditions can be found in [5] , [8] . 
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a 


M 


R a - SLP 


R a - SCP 


3.117 





1734.120 


1775.955 


3.117 


10 


1734.154 


1775.987 


3.117 


1000 


2082.802 


2100.935 


1.5 





3116.286 


31199.286 


1.5 


5 


3116.381 


3199.289 


10 





11409.157 


14909.559 


10 


100 


14414.05 


14419.963 


10 


500 


14531.694 


14994.747 


20 





166779.036 


182878.881 



Table 1 : Numerical values for the Rayleigh number for various values of the pa- 
rameters a, M obtained by spectral methods based on shifted Legendre (SLP) 
and shifted Chebyshev (SCP) polynomials. 

4 Conclusions 

In this paper we performed a bifurcation analysis and a numerical treatment 
for an electrohydrodynamic convection problem. When eigenvalue problems 
from linear stability theory are investigated only numerically spurious eigen- 
value can be encountered, especially when bifurcation sets of the characteristic 
manifold occur. In order to detect the false secular points a bifurcation study 
of the problem becomes necessary. An example of this type of problems was 
investigated in [4], e.g. for an electrohydrodynamic convection problem in the 
case of free-free boundaries the numerical methods led to the existence of false 
secular points. 

The numerical study was performed here using a Galerkin type spectral 
method which implied that the boundary conditions are satisfied by the or- 
thogonal sets of expansion functions. However, when this condition is not ful- 
filled, the tau method or the collocation method can also be applied. All these 
methods are widely used in the numerical investigation of eigenvalue problems 
governing the linear stability of motions or equilibrium of fluid in convection 
problems. From the physical point of view, the evaluations of the Rayleigh 
number R a showed an enlargement of the stability domain when the param- 
eter M is increasing, the dependence of R a of M is not however exponential. 
These evaluations were easy to compute and proved to be similar to the ones 
existing in the literature. 
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